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Introduction 


T to nredict the aerodynamic characteristics of airfoils operating at low Refolds 

numbers IR < 5 0 x 10 s ), it is necessary to accurately account for the effects of laminar (tran- 
ZSo iepar bub lies.’ In genera!, the greatest difficulty comes about when attempting 
to determine the increase in profile drag that results from the presence of separation bubbles 
Becau Te drag o„ an airfoil depends on the trailing-edge value of the momentum thickness^ 
tile increase in drag that accompanies separation bubbles is primarily due to the rap.d , pereas 
„ momentum thickness that occurs in a bubble compared to the much smaller growth tlurf 
occurs with a natural transition from a laminar to a turbulent boundary layer. £»«««; 
quence when a separation bubble is the actual mechanism of transition methods that 

transition occurs at the laminar separation point and use the boun ary- ^ underpredict 

point as the initial conditions for the turbulent boundary-layer calculations usually P 

_ Thu . to successfully determine the increase m drag on an airfoil due to separau 
bubbles! it? neiX ™tely model the development of the boundary layer through the 

bubble. 

While a number of empirically based separation bubble models have been introduced m 
the past some of which are detailed in Ref. 3, the majority of these assume that the bubo e 
development is fully predictable from upstream conditions. More recently, much more = accurate 
predictions have been made possible using viscous/inviscid interaction approaches such as tho 
bed n Refs Te. By such means, the influence of the bubble on the entire velocity d.s- 
tTutlon over the airfoil Is accounted for globally by iterating between the inviscd fiow and 
boundary-layer solutions. While not of much concern in predicting t 

teristics of a single airfoil, the amount of computational t ^^ 1 ^ = qu ; ^ed / 0 ”“^ h b t ?" n ^ ber y of 
iteration methods becomes consequential in the case of airfoil design for which 

analysis cases required can become very large. 

One way of accounting for laminar separation bubbles in airfoil design is the bubble analog 
used in the Lsign and analysis program of Eppler and Somers.’.* In this method, the dengue 
is warned about the presence of separation bubbles which might unacceptably increase ^the dra 
over that which is predicted assuming that transition occurs at laminar separation. Alth°ug 
this approach has proven very useful in designing airfoils for low Reynolds number W 1 —" . 
it would be advantageous to have predictions of section properties which more fully accoun 
for the presence of laminar separation bubbles provided this can be done without > 

increasin'* the computational time. Toward this end, a locally interactive separation bubble 
modeTiri been developed and incorporated into the Eppler and Somers program Although 

unable to account for strong interactions such as the large reduction in ““TslTusT 
caused by leading-edge bubbles, it is able to predict the increase ,n drag and the local alteration 
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of the airfoil pressure distribution that is caused by bubbles occurring in the operational range 
which is of most interest. 

To fully determine the behavior and influence of a laminar separation bubble, it is necessary 
to Juratdv preX. the shear-iayer development in the regions of the laminar separation . bubble 
indicated inF^g. 1. The formation of a bubble is initiatated at point $ shown m the figu , y 
Thetminar boundary laver separating from the airfoil surface. Using integral boundary-layer 

methods this point can be determined with the accuracy needed for airfoil design work. Once 
methods, inis po tracked and the transition from laminar to 

separated, the free shear layer development must be tracked and m separation 

turbulent flow which occurs near the point T, predicted. As shown in Fi 0 . the separa i 
bubble causes 'a plateau to form in the velocity distribution beaten the pomts — 
to laminar separation and the end of the transit, on reg, on From th s 

of the bubble encompasses a pressure recovery region which leads to the reattachment ot 

- point 1 As an additional pressure recovery "T 

of a reattachment point, the velocity distribution corresponding to the highlj non-equilibr urn, 
rllXng boundary layer downstream of reattachment -undershoots’ the mviscid distribution 
Eventually the turbulent boundary layer reaches its fully-developed state and the undershoo 
region merger smoothly from below with the inviscid velocity distr.but.on, Clearly ,t is possibly 
especially at low Reynolds numbers, that the turbulent boundary layer never reaches equilibrium 

before the trailing edge of the airfoil. 


s - LAMINAR SEPARATION POI NT 
T - [NO OF TRANSITION RE C I ON 
R - TURBULENT REATTACHWENT POINT 



Fig. 1: Sectional view of a two-dimensional short laminar separation bubble. 
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Fig. 2: Influence of a short laminar separation bubble on the velocity distri- 

bution over an airfoil. 


rwrvipw of Procedure for Mod e ler Laminar Separation Bub bles 

of the classical empirical model of Horton, modAe ^ CC ° T h reformulated in terms 

gral boundary-layer analysis employed by Eppler. Specially, th 1S method makes use 
momentum and energy integral equations, 

2) ^ a) 

ds 2 K U ds 

4¥ < 2 > 

ds U ds 

along with the appropriate closure relations . 11 

Using .he empirical separation bubble mode, “^^"TeS 

as wen," it is found WfcWW. 

only by a model which accounts for the effects of the downstream portions of the bubble on those 
Tpstrea; b particular, along with the chord Reynolds number and the upstream development 
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of tho laminar boundary layer, .he development of the bubble is found to 

P ^r:r“ ^ w « 

momentum thickness a. the laminar separation point, (*.)„ and Caster's pressure grad.en. 
parameter, 15 [(n _ l2 , rr r (AA , 1 2 A(UIU^) 


m s ) 


MUJUoo) 

A(s/c) 


It should be noted that, while (R h )s has been used in other models to account for Rasing 
transition length with increasing Reynolds number, Caster's parameter 

onlv as a criterion of when reattachment of the short bubble is not possible and a long bubble 
/ustihca^n that the value of R has a strong tnduence on the ‘--Uon pmcess an o . 
aspects of the bubble development, such as the velocity plateau region of the bubble is 

simulation presented in Ref. 16. In this research, the criterion for the 

“bursting” of short bubbles is found to correlate well with the boundary between steady and 
unsteady reattachment of a laminar boundary layer subjected to an imposed pressure ^gra _ 
The unsteadiness is due to a periodic vortex-shedding which increases m frequency and mten isiy 
L the value of Caster’s parameter increases. In view of the elliptic.* 

field, it is not too surprising that this frequency is found to correspond to the most unsta ^ 
frequency from a linear stability analysis of the laminar free shear layer. ’ 

that transition in such a Bow field would be strongly dependent on the sh.dd.ng frequMC, r « 
is characterised by the value of Caster's parameter. Because the short separ. on bubb le s 
analogous to this unsteady reattachment process, it is reasonable that m add t.on to (Kft)s. 
the transition process should depend on the downstream reattachment of the turbulent 
layer and the value of Gaster’s pressure gradient parameter. 

Startin' 1 with the inviscid velocity distribution over an airfoil, the scheme used by 
present metLd of predicting the development of laminar separation bubbles .summarized by 
the Tow diagram shown in Fig. 3. As indicated, the bubble model is evoked when laminar 
separation is predicted. At that point (R St ) S is determined and based on ^ mv lS ad vebc, y 
gradient at the laminar separation point, an initial estimate of Gas er s par , * 

From the separation point, the velocity distribution in the plateau region is prescribed The 
function which defines this distribution depends both on P and the matching o l s P 
that of the inviscid velocity distribution at the laminar separation point. Using the Prescribed 
velocity distribution over the laminar part of the bubble, the boundary- ayer eve 
accomplished using the momentum and energy integral equations Eqs. (1) and (2) ’ T flow 
appropriate closure relationships which, as will be described, were derived usmg the rever ^flow 
velocity profiles developed by Green. 17 Likewise, these profiles were used to develop the criterion 
used to predict transition in the free shear layer. This criterion is also a function of {R h )s and 
P The turbulent part of the bubble and the undershoot region are determined by prescnbin 0 
the distribution of H 22 and solving the integral boundary-layer equations in the ^ verse “° 9 
At this point, a local iteration is carried out to ensure that Horton’s reatl .achmen coition 
is satisfied for a velocity distribution in the non-equilibrium region of the bubble that merg 
smoothly downstream with the inviscid velocity distribution. Once this iteration has converged, 
a new value of P is calculated from the velocity gradient along the bubble as determined using 
the points of laminar separation and the point where the bubble recovery first crosses the mvisci 
distribution. With this new value of P, an outer iteration repeats the calculations indicated until 
the overall length of the bubble no longer changes. It is through the two iteration schemes that 
the influence of the downstream conditions in the bubble are communicated upstream. 
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Fig. 3: Flow chart of computational scheme used to model laminar separation 

bubbles in the flow over airfoils. 


Details of the T.aminar Separation Bubble Model 


In the present procedure for modelling laminar separation bubbles, the laminar closure cor- 
relations, the velocity distribution in the plateau region, transition, and the function describing 
the undershoot region depend only on {R^)s anc * ot h er features of t e u e are e- 

termined from the governing equations. The specific regions of the bubble model will now be 

described. 


Laminar Part of the Bubble 

The Eppler and Somers program uses a very reliable criterion to detect laminar separation, 
based on the value of the energy to momentum thickness shape factor, 

(tf 32 )s = 1.515095 ( 4 ) 

This value is approached from above. Since the analysis method is formulated in the direct 
mode, a small error is introduced in the prediction of the laminar separation point by the 
presence of the Goldstein singularity. This causes the distribution of H 32 to exhibit a very steep 
slope immediately upstream of separation, which is a similar behavior to that of the skm-fric ion 
coefficient. Furthermore, the direct formulation precludes any upstream influence of the bubble 
on the pressure distribution. In the present version of the bubble model this local interaction 
is neglected. Although possibly important for leading-edge bubbles, it should not have a major 
impact on the development of mid-chord bubbles. 
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Upon detection of laminar separation, the development of the separated laminar shear layer 
is calculated using the same governing equations, the momentum and energy in egra eq' “ ' 

together with closure relation, obtained from a family of reversed lamma, velocity proMe . 
Instead of implementing this boundary-layer method in the inverse mode as it is usual lly done, 
the development of a family of pressure distributions in the laminar par 
its calculation in the direct mode. 

The function used to approximate the pressure distribution in the laminar part of the bubble 
is a generalization of that developed by van Ingen and Boermans and presented m Refs. 18 and 

19, 

U f A A r A t O (5) 


where 


— = .978 + .022 exp(— 4.454£ - 2.5£ 2 ) 
Us 

3 - s$ 


{Rs 2 )s{^2)s 

This distribution, unlike the constant-pressure plateau used by Horton, allows a slight pressure 
rise after laminar separation. Using detailed pressure distributions obtained from recent win- 
tunnel tests of the NASA NLF(1)-1015 airfoil 20 in the Low-Turbulence Pressure Tunnel at the 
NASA Langley Research Center, the accuracy of Eq. (5) has been checked for several different 
conditions. It is found that, when the overall pressure gradient along the bubble decreases, e 
velocity distribution over the plateau region falls below that given by Eq. (5) while, when e 
pressure gradient steepens, it rises above. To account for this effect, Eq. (5) is rewritten as 


iU(i 

Us ^ 


DU) + DU exp(— 4.454£ - 2.5£ 2 ) 


where DU is indicative of the amount of the pressure rise accomplished by the plateau velocity 
distribution over the laminar part of the bubble. It is found to be well represented as a function 
of the Caster pressure gradient parameter, P. This functional relationship, shown m Fig. 4 
was developed by extracting corresponding values of DU and P directly from the expenmen a 
pressure distributions of the NLF(l)-1015 and the Eppler E387 airfoils. - The solid Line is a 
quadratic least-squares fit that has been included in the model, 

DU = 0.0609691 + 0.304819P -r 0.507176P 2 ( 8 ) 

It should be noted that the value of DU = 0.022 used by van Ingen and Boermans falls in the 
middle of the variation in DU shown in Fig. 4. 

Upon examination of the velocity distributions obtained using Eqs. (6)-(8), it was noticed 
that unlike in the experimental distributions, a discontinuity in the velocity gradient at the point 
of laminar separation was present. Therefore, to match the gradients at laminar separation, an 
additional variable is introduced into the velocity distribution given by Eq. (7). Specifically, 
the product \(R 63 )s (fc)s] in Eq- (6) can be treated simply as a scaling factor between the 
physical variable, s, and the universal dimensionless variable, f. Consequently, the value of this 
scaling factor is determined as that which yields a velocity distribution in the laminar part of 
the bubble whose gradient is continuous with the gradient of the inviscid velocity distribution 

at separation. 

Using the prescribed velocity distribution over the laminar part of the bubble, the develop- 
ment of the separated laminar shear layer is calculated by means of the same governing equations 
used for the rest of the boundary layer, the momentum and energy integral equations. Since 2 
and S 3 are obtained directly from the governing equations, the transition criterion is expresse 
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in terms of (tf„)r as a function of (J U,) s and P. A family of velocity profiles must likewise be 
chosen to develop the closure correlations for H 12 (^ 32 ), c/(i?32, Rs* )> ancl Gd(Hz2 » Rs , )• * 

correlations based on the reversed Falkner-Skan, or Stewartson, 22 profiles developed by Drela 
were originally included in the model. As discussed by Fitzgerald and Mueller, 23 however, recent 
measurements inside the bubble show velocity profiles that are quite different from the Stewart- 
son profiles and closer to the two-parameter profile family originally developed by Green. As 
shown in Fig. 5, the two parameters, h/b and G, are linked to the geometrical characteristics 
of the profiles, (h/b) is the ratio of the distance of the shear layer from the centerline of the 
wake to the width of the shear layer and G is the amplitude of the Coles wake function. Since 
there is slip along the centerline of such a recirculating base flow, these profiles cannot be used 
to develop a correlation for c j. In view of the characteristically small values of cj in the laminar 
part of the bubble, however, this should not be a problem. 



Fig. 4: Pressure recovery in the laminar part of the bubble as a function 

of Gaster’s pressure gradient parameter, o : NLF(l)-1015 airfoil, 
NASA LaRC, June 1987; □ : E387 airfoil. 21 


By applying the definitions for the integral thicknesses of the boundary layer and the dissi- 
pation coefficient, relationships for i?i 2 > -Hs2i C d 35 functions of h/b and G are obtaine 

for the Green profiles. The results are shown in Figs. 6 and 7 where these new two-parameter 
correlations are compared to those developed by Drela from the Stewartson profiles. The solid 
lines utilize the fitted variations of G and [h/b). As both H \2 and Hzi increase monotonically 
between separation and transition, moving to greater values of the abscissa on these plots cor- 
responds to moving downstream inside the bubble. Thus, both are similar to the Stewartson 
correlations near separation but can be quite different further downstream. It appears from 
the measurements that the back-flow, which is proportional to G, may be constant within each 
bubble although different for different bubbles. From Fig. 6, the values of shape factors actually 
measured, although different in absolute value, follow the same slope thus confirming a constant 
value of back-flow velocity. These considerations justify eliminating [h/b) between the expres- 
sions for the shape factors and expressing the closure relationships in terms of Hzii calculated 
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Fig 7: Comparison of dissipation coefficient correlations. : Ref. 5, 

from Stewartson profiles; : Ref. 23, from fitted Green profiles; 

• • • : Eq. (10). 


from the governing equations, and G, whose behavior within each bubble can be correlated to 
local flow conditions. Thus, the closure correlations are 


H 12 = 


3(1 - G) - Hz2 
(1-G)(1 — 2G) 


Rs^Cd = 


7T 2 G 3 

2 


3 _ (4 - 5G)(1 - G) - (2 - 3G)i?32 

2^ 4(1 — G) — 2Hzz 


( 9 ) 

( 10 ) 


The unknown quantities in the boundary-layer method for the laminar part of the bubble 
have been reduced to G and c/. Physically grounded assumptions can be made about the 
dependence of these variables on local flow conditions. For instance, it seems reasonable to 
expect that, as the pressure gradient along the bubble increases, so does the strength of the 
recirculation and, therefore, — Cf. In the present version of the model, however, G is related 
to P such that (S 2 ) s plays an important part in determining the ratio of reverse velocity to 
forward boundary-layer edge velocity. The value of cj is held constant between separation and 
transition. 


Transition 


The criterion for predicting transition in the free shear layer presently used was arrived at 
only after establishing that the addition of a second parameter in the laminar reversed veloc- 
ity profiles still did not provide the bubble model with enough generality. The Green profiles 
correlations have been retained, however, as the present form of the transition criterion is be- 
lieved to embody both the stability characteristics of free shear layers of varying velocity ratios 
(varying amounts of back-flow) as well as the perhaps more important input from the unsteady 
reattaching turbulent shear layer. 
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6: Comparison of shape factor correlations. — : Ref- a, R° m 

Stewartson profiles; : Ref. 23, from fitted Green profiles, 

... : Eq. (9): o : measured, Ref. 23. 



In the present separation bubble model, transition is 
shape factor, (# 32 ) 7 , is achieved. As with other features 
dependent on both (i?5 3 )s and P. This relationship is 


predicted when a given value of the 
of the bubble, (^ 32)7 is taken to be 


[H Z i)r 


875 - (iZfe)s 

35000000P 2 + 1800000P -r 40000 


+ 1.515095 


( 11 ) 


To explain the transition criterion, it is plotted together with the shear layer development on e 
same plot that Eppler uses to describe the boundary-layer development. Since #32 and R 6 , are 
calculated at each point along the boundary layer, by connecting subsequent (P 32 , PfJ-pairs 
on a plot whose axes correspond to these two variables, the boundary-layer developmen rom 
the stagnation point to the trailing edge can be described in a very concise way. ig* s ows 
one such boundary-layer development together with all the transition and separation criteria. 


LOG 
R 0 2 



Fig. 8: Modified Eppler boundary-layer development plot showing the 

boundary-layer development inside the bubble. 


Turbulent Part of the Bubble and Reattachment 

An analysis recently completed by R. Eppler that establishes the greatest possible pressure 
recovery for an attached turbulent boundary layer provided inspiration for developing a new 
approach to the treatment of the turbulent part of the bubble. Given that the value as well as 
the slope of the H 27 distribution is always known at the transition and reattachment points, 
a general function has been developed which allows the solution of the turbulent part of bhe 
bubble in the inverse mode. The distribution of H 32 I s specified as 
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H Z2 {o) = 


(g.«)T+(Hii)» + Uhill- _ Uhil *. sin 0 - a - 1 


(F 32 )je + ^ 


1 + sin 


(^+•01 


o > 1 


where 


cr = 


5 - sy 
ti 


The corresponding distributions of C/(.) and S, are obtained from 

dU Hz 7 c f ~~ Cd + ^ 2-^32 


ds 62^32(^12 “ l) 


■U 


d 6 2 
ds 




( 12 ) 


(13) 


(14) 


(15) 


H n - 1 ' (-H 12 - 1)#32 

where//' denotes the derivative with respect to s of Eq. (12). The closure relations for 

Teparated a'l attached turbulent boundary layers developed by Drela’ are used, enforcing 
continuity in i?32 transition. 

In the expression defining the reat.achment process, Eq. (12) there are two 1 free constants 
. „ d A such that two additional conditions are necessary to achieve a unique solution. One 
condition' is that the undershoot merge smoothly with the inviscid velocity distribution. Simul- 
taneously, the reattachment condition derived by Horton' must be satisfied, 

(16) 




Cd 1 

Udx 

z 

.# 32(^12 “ 1 ) . 


The accuracy of the self-consistent reattachment process described above depends ; on the 
choice of the function that represents the input if 3! distribution ind ‘he value f (Has)^ 
As the value of Ha, at the end of the undershoot, which depends on the amplitude A, does 
not vary significantly with Reynolds number or pressure gradient, the function employed is 
Wiev^ t. provide a satisfactory approximation to the actual houndary-la^r d^opm« m 

the reattachment region as it is affected by the local interaction. Therefore the only pa ame.er 

requiring further calibration is (Airly. The transition criterion, Eq. (11 , can be J 
matching predicted and measured pressure distributions. Havmg fixed (Ass)r, the correct 
dependence of G on P and, possibly, on <R h )s, can be determined by matching predicted and 

measured transition lengths. 

Prediction 1 ^ Rased on the Present Bub ble Mode] 

At the present time, all of the components of the bubble model are in place with enough 
flexibility to capture the physics for a wide range of conditions. It remains, however to determine 
the exact forms of the correlation functions that will “tune” the model to achieve the gene^bty 
and robustness that this wide range requires. With these qualifications in mind, the pressur 
distribution for the Eppler E387 airfoil at an angle of attack of 1.5 degrees as predicted isi^g 
the Eppler and Somers program with the bubble model included 15 s ° wn m 6 ‘ . , 

included in the figure is the distribution obtained with the Drela and Giles program 
critical amplification factor of 12, along with that of McGhee et al.- .both at 2 degrees, 
this case, the predicted pressure distributions agree reasonably well with one another and 
the experimental results. As expected, because this airfoil is not highly aft-loaded, the influence 
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Fi.. 9: Velocity distribution and upper-surface boundary-layer development 

predicted by present method compared to those of the Drela and Giles 

program.* Eppler E387 airfoil, c, = .55. — : present method, 

" - 1 — 91 




of viscous/inviscid interactions is not strong. Consequently, outside of the bubble regi 
inviscid prediction is comparable to that of the fully interactive method and, within the bubble 
region, the local interactive calculation of the present model produces results which are very 
close to those of the fully interactive method. 

A comparison of the development of the boundary-layer properties for the E387 airfoil, H u , 
S 2 , H * 2 , and C/ , calculated using the present method with those obtained using the re a an 
Giles program are also presented in Fig. 9. As with the pressure distnl ^ tl ° n ^ t , h % a ^ ee ” e e n r 
between the predictions obtained using the present bubble model with those of the fully inter- 
active method is good. The apparent difference in bubble length may be due to the different 
treatment of the transition region between the two methods. As in most transition prediction 
methods based on linear stability, the critical amplification factor used m t e re a an 1 e 
program marks the beginning of the transition region while, in the present model the transi- 
tion point used corresponds to the end of this region and the beginning of the full turbulent 
calculations. The difference in the distribution of H Z2 can be attributed to the 1 erence m 
pressure distributions in the laminar part of the bubble, shown in Fig. 9 , as predicted by the 
two methods. The steeper pressure gradient predicted by the present model leads to a larger 

growth in S 2 and, therefore, a smaller H 33 - 
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The aerodynamic charncCerisucs for the E387 ai« 1 

it. In fa , , w be noted that because the original program assumes transition 

“The toin« separation point, the width between the original transition development and that 
Obtained with the modified program is the length of the laminar part of the bubble. 


387 


9.0656 


T. - boundary layer transition 
S. = boundary layer seporotion 

U. ■= upper surfoce 
L. - lower surface 



Because of the strong viscous interaction due to the large amount of 
the NLF(1)-1015 airfoil represents a much more difficult case for the presen me w 

the E387 airfoil. The pressure distribution for this airfoil at an angle 
the chordline of -3.0 degrees is shown in Fig. 11. The greater impact viscous 
in this case compared to the last is apparent. While it should be emphasised again ‘hat the 
present model has yet to be fully calibrated, it is seen in this case that the overall method doe 
seem to model the physical phenomenon sufficiently well to capture the two different bubble 
developments that occur on this airfoil at this angle of attack. 

Concluding Remarks 

While its development is still in progress, the laminar separation bubble model described 
seems to properly predict the behavior of the various parts of the bubble and possess e * “ 
ibility required to work over a wide range of cases. For some cases, such as those m wbch he 
viscous interactions greatly modify the inviscid pressure distribution, high accuracy in pr 
ing the drag and the influence of the bubble clearly requires a fully interactive method. For 
many other cases, however, the present local interactive bubble model should provide a C ° ? ' 
tationally efficient method which is well suited for the aerodynamic analysis which is required 

during the process of airfoil design. 


13 




0.0 0.5' x/c 1*0 

RE ( 6 ) = 0.500 AIPHA= -4.00 


Fig 11: Velocity distribution and boundary-layer development predicted by 

present method for NLF(l)-1015 airfoil, c* = .27. o : experimental, 
NASA LaRC, June 1987. 
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